Task 2¶

Integrantes¶

  • Sergio Orellana 221122
  • Rodrigo Mansilla 22611
  • Ricardo Chuy 221007

El objetivo de esta parte es implementar el pipeline de alineación sin depender de cajas negras. Por ello considere que no deben de usar cv2.findHonography o cv2.RANSAC. Además para este laboratorio necesitara crear su propio dataset, por ello tome 3 fotografías propias de una escena planar (i.e. una pancarta en una pared, un cuadro, o una fachada de edificio lejana) con ángulos y perspectivas drásticamente diferentes.

In [2]:
import cv2

import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np

img1 = mpimg.imread('imgs/abajo.jpg')
img2 = mpimg.imread('imgs/frente.jpg')
img3 = mpimg.imread('imgs/lado.jpg')
# por algun motivo la imagegn 3 sale horizontalmente, entonces le di vuelta jej
img3 = np.rot90(img3, k=-1)

fig, axes = plt.subplots(1, 3, figsize=(15, 5))  # 1 fila, 3 columnas

axes[0].imshow(img1)
axes[0].set_title('Vista de abajo')
axes[0].axis('off')  

axes[1].imshow(img2)
axes[1].set_title('Vista de frente')
axes[1].axis('off')

axes[2].imshow(img3)
axes[2].set_title('Vista de lado')
axes[2].axis('off')

# plt.tight_layout()  # Ajustar espaciado
plt.show()
No description has been provided for this image

1. Detección y Macheo¶

  • a. Utilice SIFT u ORB (permitido usar OpenCV aquí) para detectar puntos de interés y descriptores.
  • b. Realice un emparejamiento de fuerza bruta (Brute-Force Matcher).
  • c. Requisito: Visualice los matches antes de filtrar. Debe verse una cantidad considerable de ruido/errores

Parte a¶

In [3]:
PATH_REF = 'imgs/frente.jpg'   
PATH_OBJ = 'imgs/abajo.jpg'     

print(f"Referencia: {PATH_REF}\n Objetivo:   {PATH_OBJ}")

# 
img_ref_bgr = cv2.imread(PATH_REF)
img_obj_bgr = cv2.imread(PATH_OBJ)

# Validación simple
if img_ref_bgr is None or img_obj_bgr is None:
    raise ValueError("No se encontro la img")

img_ref_rgb = cv2.cvtColor(img_ref_bgr, cv2.COLOR_BGR2RGB)
img_obj_rgb = cv2.cvtColor(img_obj_bgr, cv2.COLOR_BGR2RGB)

img_ref_gray = cv2.cvtColor(img_ref_bgr, cv2.COLOR_BGR2GRAY)
img_obj_gray = cv2.cvtColor(img_obj_bgr, cv2.COLOR_BGR2GRAY)
sift = cv2.SIFT_create()

# los keypoints y descriptores
kp_ref, des_ref = sift.detectAndCompute(img_ref_gray, None)
kp_obj, des_obj = sift.detectAndCompute(img_obj_gray, None)

print(f"Keypoints en Referencia: {len(kp_ref)}")
print(f"Keypoints en Objetivo:   {len(kp_obj)}")

Referencia: imgs/frente.jpg
 Objetivo:   imgs/abajo.jpg
Keypoints en Referencia: 28311
Keypoints en Objetivo:   10239

Parte b¶

In [4]:
# Usamos Fuerza Bruta con distancia Euclidiana, el k=2 igual que en el lab pasado significa
# los 2 mejores matches
matcher = cv2.BFMatcher(cv2.NORM_L2, crossCheck=False)
matches_crudos = matcher.knnMatch(des_ref, des_obj, k=2)

print(f"Matches Totales: {len(matches_crudos)}")

Matches Totales: 28311

Parte c¶

In [5]:
matches_sin_filtro = [m for m, n in matches_crudos]

img_ruido = cv2.drawMatches(
    img_ref_rgb, kp_ref, 
    img_obj_rgb, kp_obj, 
    matches_sin_filtro, # Pasamos la muestra ruidosa
    None, 
    flags=2
)

plt.figure(figsize=(15, 6))
plt.title("Matches Crudos sin filtro")
plt.imshow(img_ruido)
plt.axis('off')
plt.show()
No description has been provided for this image
In [6]:
# esto es algo extra que quisimos poner para verificar el match con filtro como el lab pasado

ratio = 0.75
matches_buenos = []

for m, n in matches_crudos:
    if m.distance < ratio * n.distance:
        matches_buenos.append(m)

print(f"Matches Totales: {len(matches_crudos)}")
print(f"Matches Buenos: {len(matches_buenos)}")

# visualizacion de los matches
img_matches = cv2.drawMatches(img_ref_rgb, kp_ref, img_obj_rgb, kp_obj, matches_buenos, None, flags=2)
plt.figure(figsize=(15, 5))
plt.title("Matches Buenos")
plt.imshow(img_matches)
plt.axis('off')
plt.show()

Matches Totales: 28311
Matches Buenos: 1102
No description has been provided for this image

2. Algoritmo DLT¶

a. Escriba una función calcular_homografia_dlt(puntos_src, puntos_dst) que reciba exactamente 4 pares de puntos.

b. Debe construir la matriz A de tamaño 8×9.

c. Debe resolver el sistema usando SVD (numpy.linalg.svd).

d. Nota: Debe normalizar los puntos antes de entrar al DLT (restar la media y dividir por la desviación estándar) para estabilidad numérica, y des-normalizar la matriz H resultante.

Función de Homografía¶

In [7]:
def normalizar_puntos(puntos, eps=1e-12):
    puntos = np.asarray(puntos, dtype=np.float64)
    if puntos.ndim != 2 or puntos.shape[1] != 2:
        raise ValueError("puntos debe tener forma (N,2)")

    mu = np.mean(puntos, axis=0)
    sigma = np.std(puntos, axis=0)
    sigma[sigma < eps] = 1.0  

    T = np.array([
        [1.0/sigma[0], 0.0, -mu[0]/sigma[0]],
        [0.0, 1.0/sigma[1], -mu[1]/sigma[1]],
        [0.0, 0.0, 1.0]
    ], dtype=np.float64)

    pts_h = np.hstack([puntos, np.ones((puntos.shape[0], 1))])
    pts_n = (T @ pts_h.T).T
    return pts_n[:, :2], T

def calcular_homografia_dlt(puntos_src, puntos_dst, eps=1e-12):
    """
    DLT con N>=4 pares de puntos (para refinamiento con inliers).
    Construye la matriz A de tamaño 2N×9 y resuelve con SVD.
    """
    src = np.asarray(puntos_src, dtype=np.float64)
    dst = np.asarray(puntos_dst, dtype=np.float64)
    
    if src.shape[0] < 4 or dst.shape[0] < 4:
        raise ValueError("Se requieren al menos 4 pares de puntos.")
    if src.shape[0] != dst.shape[0] or src.shape[1] != 2 or dst.shape[1] != 2:
        raise ValueError("src y dst deben tener forma (N,2) con el mismo N.")

    src_n, T_src = normalizar_puntos(src)
    dst_n, T_dst = normalizar_puntos(dst)

    # Matriz A de 2N×9
    A = []
    for (x, y), (xp, yp) in zip(src_n, dst_n):
        A.append([-x, -y, -1,  0,  0,  0, x*xp, y*xp, xp])
        A.append([ 0,  0,  0, -x, -y, -1, x*yp, y*yp, yp])

    A = np.asarray(A, dtype=np.float64)  # (2N, 9)
    _, _, Vt = np.linalg.svd(A)
    Hn = Vt[-1].reshape(3, 3)

    H = np.linalg.inv(T_dst) @ Hn @ T_src
    if abs(H[2, 2]) < eps:
        return H
    return H / H[2, 2]

3. RANSAC Manual¶

a. Implemente la función ransac_homografia(matches, umbral, prob_exito).

b. Cálculo de N: Su código debe calcular dinámicamente el número de iteraciones N basado en la fórmula de probabilidad vista en clase. No "hardcodee" el número 1000.

c. Bucle: i. Seleccione 4 matches aleatorios. ii. Llame a su función DLT. iii. Proyecte todos los puntos fuente usando H_test. iv. Calcule el error de reproyección (distancia Euclidiana) y cuente los inliers.

a. Refinamiento: Una vez encontrado el mejor conjunto de inliers, recalcule 𝐻 final usando todos los inliers (no solo los 4 iniciales) mediante SVD

In [8]:
def ransac_homografia(matches, kp_src, kp_dst, umbral, prob_exito=0.99):
    """
    RANSAC manual para estimar homografía OBJ → REF.
    
    matches fueron creados con knnMatch(des_ref, des_obj),
    por lo que queryIdx → kp_ref y trainIdx → kp_obj.
    
    Para estimar H: OBJ → REF directamente, se usa:
      src_pts = puntos en OBJ (trainIdx)
      dst_pts = puntos en REF (queryIdx)
    """
    obj_pts = np.array([kp_dst[m.trainIdx].pt for m in matches])   
    ref_pts = np.array([kp_src[m.queryIdx].pt for m in matches])   
    num_matches = len(matches)

    s = 4               # puntos por muestra
    w = 0.5             # estimación inicial de ratio de inliers
    N = int(np.ceil(np.log(1 - prob_exito) / np.log(1 - w**s))) if w < 1 else 1

    mejor_H = None
    mejor_inliers = []
    iteracion = 0

    while iteracion < N:
        if num_matches < 4:
            break
        idx = np.random.choice(num_matches, size=4, replace=False)
        
        src_4 = obj_pts[idx]   
        dst_4 = ref_pts[idx]   

        try:
            H_test = calcular_homografia_dlt(src_4, dst_4)  # DLT con exactamente 4 pts (A=8×9)
        except:
            iteracion += 1
            continue

        pts_h = np.hstack([obj_pts, np.ones((num_matches, 1))])  
        proj_h = (H_test @ pts_h.T).T
        proj = proj_h[:, :2] / proj_h[:, 2:3]

        # distancia Euclidiana
        errores = np.linalg.norm(proj - ref_pts, axis=1)
        inliers = np.where(errores < umbral)[0]

        if len(inliers) > len(mejor_inliers):
            mejor_inliers = inliers
            mejor_H = H_test

            # Actualizar w y recalcular N dinámicamente
            w = len(inliers) / num_matches
            if w == 1:
                N = 1
            else:
                N = int(np.ceil(np.log(1 - prob_exito) / np.log(1 - w**s)))

        iteracion += 1

    # Refinamiento: recalcular H con TODOS los inliers
    if len(mejor_inliers) >= 4:
        H_final = calcular_homografia_dlt(
            obj_pts[mejor_inliers],
            ref_pts[mejor_inliers]
        )
    else:
        H_final = mejor_H

    print(f"Inliers: {len(mejor_inliers)} / {num_matches}")
    print(f"Iteraciones: {iteracion}")
    return H_final, mejor_inliers

Salida¶

In [9]:
# ===== Celda: Salida =====

H_final, inliers = ransac_homografia(
    matches_buenos, kp_ref, kp_obj, 
    umbral=5.0, prob_exito=0.99
)

print("Homografía encontrada (OBJ → REF):")
print(H_final)

# H_final ya es OBJ → REF, no necesitamos inv()
h_ref, w_ref = img_ref_rgb.shape[:2]
h_obj, w_obj = img_obj_rgb.shape[:2]

# Proyectar las esquinas de OBJ al espacio de REF
esquinas_obj = np.array([[0,0],[w_obj,0],[w_obj,h_obj],[0,h_obj]], dtype=np.float64)
esquinas_h = np.hstack([esquinas_obj, np.ones((4, 1))])
esquinas_proj = (H_final @ esquinas_h.T).T
esquinas_proj = esquinas_proj[:, :2] / esquinas_proj[:, 2:3]

todas_esquinas = np.vstack([
    esquinas_proj, 
    [[0,0],[w_ref,0],[w_ref,h_ref],[0,h_ref]]
])

x_min = int(np.floor(todas_esquinas[:, 0].min()))
y_min = int(np.floor(todas_esquinas[:, 1].min()))
x_max = int(np.ceil(todas_esquinas[:, 0].max()))
y_max = int(np.ceil(todas_esquinas[:, 1].max()))

T_traslacion = np.array([[1,0,-x_min],[0,1,-y_min],[0,0,1]], dtype=np.float64)
ancho_canvas = x_max - x_min
alto_canvas = y_max - y_min

img_warp_obj = cv2.warpPerspective(
    img_obj_rgb, 
    T_traslacion @ H_final, 
    (ancho_canvas, alto_canvas)
)

mask_ones = np.ones((h_obj, w_obj), dtype=np.uint8) * 255
mask_warp_obj = cv2.warpPerspective(
    mask_ones, 
    T_traslacion @ H_final, 
    (ancho_canvas, alto_canvas)
)
mask_obj = (mask_warp_obj > 0).astype(np.float64)

img_warp_ref = np.zeros((alto_canvas, ancho_canvas, 3), dtype=np.uint8)
img_warp_ref[-y_min:-y_min+h_ref, -x_min:-x_min+w_ref] = img_ref_rgb

mask_ref = np.zeros((alto_canvas, ancho_canvas), dtype=np.float64)
mask_ref[-y_min:-y_min+h_ref, -x_min:-x_min+w_ref] = 1.0

mask_obj_3 = np.stack([mask_obj]*3, axis=2)
mask_ref_3 = np.stack([mask_ref]*3, axis=2)
total = mask_obj_3 + mask_ref_3
total[total == 0] = 1.0

panorama = (img_warp_obj.astype(np.float64)*mask_obj_3 + 
            img_warp_ref.astype(np.float64)*mask_ref_3) / total
panorama = np.clip(panorama, 0, 255).astype(np.uint8)

# Visualización
fig, axes = plt.subplots(2, 2, figsize=(18, 14))
axes[0,0].imshow(img_obj_rgb);  axes[0,0].set_title("Original: Vista de abajo"); axes[0,0].axis('off')
axes[0,1].imshow(img_ref_rgb);  axes[0,1].set_title("Referencia: Vista frontal"); axes[0,1].axis('off')
axes[1,0].imshow(img_warp_obj); axes[1,0].set_title("Imagen de abajo warpeada"); axes[1,0].axis('off')
axes[1,1].imshow(panorama);     axes[1,1].set_title("Panorama"); axes[1,1].axis('off')
plt.tight_layout()
plt.show()

Inliers: 739 / 1102
Iteraciones: 21
Homografía encontrada (OBJ → REF):
[[ 2.06100658e+00  7.35563316e-01 -1.17781016e+03]
 [-3.35760489e-02  3.79289283e+00 -9.98153206e+00]
 [-4.34022351e-05  7.20384772e-04  1.00000000e+00]]
No description has been provided for this image

Task 3¶

Integrantes¶

  • Sergio Orellana 221122
  • Rodrigo Mansilla 22611
  • Ricardo Chuy 221007

1. Evaluación del umbral de RANSAC¶

a) Gráfica: Umbral vs Número de Inliers¶

In [10]:
import cv2
import numpy as np
import matplotlib.pyplot as plt

img_ref_bgr = cv2.imread('imgs/frente.jpg')
img_obj_bgr = cv2.imread('imgs/abajo.jpg')

img_ref_rgb = cv2.cvtColor(img_ref_bgr, cv2.COLOR_BGR2RGB)
img_obj_rgb = cv2.cvtColor(img_obj_bgr, cv2.COLOR_BGR2RGB)

img_ref_gray = cv2.cvtColor(img_ref_bgr, cv2.COLOR_BGR2GRAY)
img_obj_gray = cv2.cvtColor(img_obj_bgr, cv2.COLOR_BGR2GRAY)

sift = cv2.SIFT_create()
kp_ref, des_ref = sift.detectAndCompute(img_ref_gray, None)
kp_obj, des_obj = sift.detectAndCompute(img_obj_gray, None)

matcher = cv2.BFMatcher(cv2.NORM_L2, crossCheck=False)
matches_crudos = matcher.knnMatch(des_ref, des_obj, k=2)

ratio = 0.75
matches_buenos = [m for m, n in matches_crudos if m.distance < ratio * n.distance]

umbrales = [0.5, 1, 2, 3, 5, 8, 10, 15, 20, 30, 50]
num_inliers_list = []
num_iteraciones_list = []

np.random.seed(42)  

for umbral in umbrales:
    H, inliers = ransac_homografia(
        matches_buenos, kp_ref, kp_obj,
        umbral=umbral, prob_exito=0.99
    )
    num_inliers_list.append(len(inliers))
    print(f"  Umbral: {umbral:5.1f} px  →  Inliers: {len(inliers)}")
    print()

fig, ax = plt.subplots(figsize=(10, 6))

ax.plot(umbrales, num_inliers_list, 'bo-', linewidth=2, markersize=8)

for u, n in zip(umbrales, num_inliers_list):
    ax.annotate(f'{n}', (u, n), textcoords="offset points", 
                xytext=(0, 12), ha='center', fontsize=9)

ax.set_xlabel('Umbral de error (píxeles)', fontsize=12)
ax.set_ylabel('Número de Inliers', fontsize=12)
ax.set_title('RANSAC: Umbral vs Número de Inliers', fontsize=14)
ax.set_xscale('log')
ax.set_xticks(umbrales)
ax.set_xticklabels([str(u) for u in umbrales])
ax.grid(True, alpha=0.3)

ax.axhline(y=len(matches_buenos), color='r', linestyle='--', alpha=0.5, 
           label=f'Total matches buenos: {len(matches_buenos)}')
ax.legend(fontsize=11)

plt.tight_layout()
plt.show()

C:\Users\rodri\AppData\Local\Temp\ipykernel_27840\1035500192.py:40: RuntimeWarning: divide by zero encountered in divide
  proj = proj_h[:, :2] / proj_h[:, 2:3]
C:\Users\rodri\AppData\Local\Temp\ipykernel_27840\1035500192.py:40: RuntimeWarning: invalid value encountered in divide
  proj = proj_h[:, :2] / proj_h[:, 2:3]

Inliers: 148 / 1102
Iteraciones: 14154
  Umbral:   0.5 px  →  Inliers: 148

Inliers: 349 / 1102
Iteraciones: 456
  Umbral:   1.0 px  →  Inliers: 349

Inliers: 544 / 1102
Iteraciones: 76
  Umbral:   2.0 px  →  Inliers: 544

Inliers: 667 / 1102
Iteraciones: 45
  Umbral:   3.0 px  →  Inliers: 667

Inliers: 671 / 1102
Iteraciones: 32
  Umbral:   5.0 px  →  Inliers: 671

Inliers: 773 / 1102
Iteraciones: 17
  Umbral:   8.0 px  →  Inliers: 773

Inliers: 842 / 1102
Iteraciones: 12
  Umbral:  10.0 px  →  Inliers: 842

Inliers: 853 / 1102
Iteraciones: 11
  Umbral:  15.0 px  →  Inliers: 853

Inliers: 854 / 1102
Iteraciones: 11
  Umbral:  20.0 px  →  Inliers: 854

Inliers: 863 / 1102
Iteraciones: 10
  Umbral:  30.0 px  →  Inliers: 863

Inliers: 867 / 1102
Iteraciones: 10
  Umbral:  50.0 px  →  Inliers: 867
No description has been provided for this image
In [ ]:

b) Discusión: Riesgos del umbral¶

Umbral demasiado estricto (ej. 0.5px):

Respuesta:

Umbral demasiado laxo (ej. 50px):

Respuesta:

2. Caso práctico: Drones e imágenes térmicas de paneles solares¶

a) Contexto del problema¶

El drone vuela a 50 metros de altura. El terreno no es perfectamente plano (hay colinas suaves), pero los paneles sí son planos.

b) Pregunta A: ¿Es válido usar una Homografía global para unir todo el mapa del terreno?¶

Respuesta:

c) Pregunta B: Estrategia para reducir tiempo de RANSAC¶

Contexto: RANSAC tarda 3 segundos por frame en Raspberry Pi. El 90% de los matches iniciales son outliers debido al pasto y árboles repetitivos.

Propuesta concreta:

Respuesta: